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Abstract We present results of non-linear, 2D, numerical simulations of magneto- 
acoustic wave propagation in the photosphere and chromosphere of small-scale flux 
tubes with internal structure. Waves with realistic periods of three to five minutes 
are studied, after applying horizontal and vertical oscillatory perturbations to the 
equilibrium model. Spurious reflections of shock waves from the upper boundary are 
minimized thanks to a special boundary condition. This has allowed us to increase 
the duration of the simulations and to make it long enough to perform a statistical 
analysis of oscillations. The simulations show that deep horizontal motions of the 
flux tube generate a slow (magnetic) mode and a surface mode. These modes are 
efficiently transformed into a slow (acoustic) mode in the va < cs atmosphere. 
The slow (acoustic) mode propagates vertically along the field lines, forms shocks 
and remains always within the flux tube. It might deposit effectively the energy of 
the driver into the chromosphere. When the driver oscillates with a high frequency, 
above the cut-off, non-linear wave propagation occurs with the same dominant driver 
period at all heights. At low frequencies, below the cut-off, the dominant period of 
oscillations changes with height from that of the driver in the photosphere to its first 
harmonic (half period) in the chromosphere. Depending on the period and on the 
type of the driver, different shock patterns are observed. 

Keywords: Photosphere, Chromosphere, Magnet ohydrodynamics. Oscillations, Mag-| 
netic fields 



1. Introduction 

Solar magnetic structures show a continuous distribution of fluxes and sizes. In active 
regions, the field strength concentrated in large-scale structures (sunspots and pores) 
can be as large as 2 — 4 kG. In plage and network areas, the average flux decreases 
to several hundreds Gauss. However, individual flux tubes in these areas can have 
intrinsic magnetic field strength similar to that of solar pores, i.e. 1 — 2 kG. The 
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decrease of the average flux is caused mainly by the decrease of the size of the 
magnetic features, i.e. their filling factor. 

The observed photospheric brightness of the magnetic structures is in close re- 
lationship with their flux and size (see for example Berger, Rouppe van der Voort, 
and Lofdahl, 2007). Large-scale features like sunspots and pores appear dark in 
photospheric observations. Plage and facular flux tubes appear as bright features 
(see for example Figure 7 in Berger, Rouppe van der Voort, and Lofdahl, 2007). 

Due to their larger size, sunspots show horizontal gradients of the magnetic field 
strength and inclination. Instead, small-scale flux tubes possess a more homogeneous 
field with a rather sharp transition between the magnetized and non-magnetized sur- 
roundings. From the point of view of wave propagation, it is important to realize that 
the temperature at a given geometrical height is different in the different structures. 
Temperature and magnetic field define the characteristic wave propagation speeds, 
i.e. sound speed, c^, and Alfven speed, va, as well as the acoustic cut-off frequency, 
cue 5 and the height of the transformation layer where cs is equal to va and different 
wave modes can interact. These parameters can be rather different in large-scale 
dark structures and small-scale bright structures. Thus, it is not surprising that 
the observed properties of waves in sunspots, network, and plage areas are rather 
different. 

1.1. Sunspots and Pores 

Several decades of studies of sunspot oscillations can be summarized as follows. 
Waves observed at different layers of the atmosphere of sunspots umbrae and penum- 
brae seem to be the manifestation of the same phenomenon (see e.g. Malt by et al.^ 

1999, 2001, Brynildsen et a/., 2000, 2002, Christopoulou, Georgakilas, and Koutchmy, 

2000, 2001, Rouppe van der Voort et al, 2003, Tziotziou et a/., 2006). In the pho- 
tosphere, five-minute oscillations are observed in sunspot umbrae and penumbrae 
and in solar pores (Bogdan and Judge, 2006). Wave propagation is linear with 
amplitudes of several hundreds of m s""*^ and is parallel to the magnetic field lines 
(see e.g. Gurman and Leibacher, 1984, Lites, 1984, Collados et a/., 2001, Centeno, 
Collados, and Trujillo Bueno, 2006, Bloomfield, Lagg, and Solanki, 2007). Magnetic 
field oscillations of a few Gauss are detected in the deeper layers and attempts have 
been made to interpret these oscillations in terms of fast and slow MHD waves {e.g. 
Riiedi et a/., 1998, Lites et al, 1998, Bellot Rubio et a/., 2000, Khomenko, Collados, 
and Bellot Rubio, 2003). 

In the chromosphere, the dominant period of oscillations decreases and shock 
waves are observed with a three-minute periodicity in sunspot umbrae and in pores 
(Lites, 1984, 1986, 1988, Socas-Navarro, Trujillo Bueno, and Ruiz Cobo, 2000, Cen- 
teno, Collados, and Trujillo Bueno, 2006 Tziotziou et a/., 2007). The amplitudes of 
the observed shocks are in the range of 5 — 15 km s""*^, depending on the size of the 
structure. In pores, shocks are weaker (Centeno, Collados, and Trujillo Bueno, 2006). 
Umbral flashes also show a three-minute periodicity. In contrast, in the penumbra, 
the shock behavior of waves with five-minute periodicity is detected for running 
penumbral waves even at chromospheric heights (Tziotziou et a/., 2006). Waves ob- 
served in the transition region and the corona conserve these properties (De Moortel 
et al, 2002, Marsh and Walsh, 2005, 2006, Maltby et al, 1999, 2001, Brynildsen 
et a/., 2000, 2002, Christopoulou, Georgakilas, and Koutchmy, 2000, 2001). Wave 
propagation from the photospheric pulse into the chromosphere and higher layers is 
along the magnetic field lines for both three and five minutes perturbations. 
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1.2. Plages, Facular regions, and Network 

Bright magnetic field structures like plages, facular regions, and network manifest 
five-minute linear oscillations in the photosphere and the oscillations maintain the 
five minute periodicity in the chromosphere and larger heights (Lites, Rutten, and 
Kalkofen, 1993, Krijger et a/., 2001, De Pontieu, Erdelyi, and de Wijn, 2003, Centeno, 
Collados, and Trujillo Bueno, 2006, Bloomfield et a/., 2006, Vecchio et a/., 2007). 
Shock occurrence is much lower in these quiet Sun regions with enhanced flux (Krijger 
et a/., 2001, Centeno, Collados, and Trujillo Bueno, 2006). The velocity oscillation 
amplitudes in the chromosphere do not reach values as large as in sunspots and stay 
within 2 — 4 km s~^. Waves at photospheric and chromospheric layers seem to be cor- 
related and vertical propagation along the vertical magnetic fields seem to dominate 
(Krijger et a/., 2001, Centeno, Collados, and Trujillo Bueno, 2006, Vecchio et a/., 
2007). Curiously, there are areas with enhanced three-minute power surrounding the 
strong magnetic field concentrations, possibly related to flux tube canopies (aureoles, 
see Braun and Lindsey, 1999, Thomas and Stanchfield, 2000, Krijger et a/., 2001, 
Judge, Tarbell, and Wilhelm, 2001). In the higher coronal layers, loops connecting 
sunspots oscillate with three minutes and the loops connecting network and plage 
regions oscillate with five minute periods (De Moortel et a/., 2002). 

Thus, the wave properties are distinctly different in large-scale and small-scale 
magnetic features. The difference in oscillation properties is present through the 
whole solar atmosphere and should be a consequence of the thermal and magnetic 
properties of these structures. 

1.3. Theoretical Models 

A number of theoretical mechanisms have been proposed to explain the oscillation 
spectra at different heights in magnetic structures and to interpret the observed 
oscillations in terms of MHD waves. Zhugzhda and Locans (1981), Gurman and 
Leibacher (1984), and Zhugzhda (2007) argue that the observed spectrum of sunspot 
umbral oscillations in the chromosphere is due to the temperature gradients of the 
atmosphere acting as an interference filter for linear three-minute period acoustic 
waves. On the other hand, it has been demonstrated by Fleck and Schmitz (1991) 
that the change of the period with height from five to the three minutes is a basic 
phenomena occurring even in an isothermal atmosphere for linear waves due to 
the resonant excitation at the atmospheric cut-off frequency. In the case of the 
solar atmosphere, the temperature minimum gives rise to the three-minute cut-off 
period. Later, the same authors found that, due to non-linear shock interaction and 
overtaking, the frequency of acoustic oscillations also changes with height (Fleck 
and Schmitz, 1993). The response of the solar atmosphere to an adiabatic shock 
wave leads to a natural appearance of the three-minute peak in the oscillation power 
spectra, in the case the underlying photosphere has a five-minute periodicity. A 
third mechanism was shown to explain the three-minute periodicity of calcium bright 
grains in the quiet Sun: with a spectrum of periods in the photosphere with a peak at 
five minutes, this distribution changes to one with a maximum at the acoustic cut-off 
because at longer periods the energy falls off exponentially with height (Carlsson and 
Stein, 1997). These mechanisms may explain the shift with height of the period of 
oscillations from five to three minutes observed in sunspots. 
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The question remains: why is there not such a shift in smah-scale structures of 
plages and network regions? De Pontieu, Erdelyi, and Stewart (2004) argue that 
the inclination of the magnetic field may play an important role. If (acoustic or 
slow MHD) waves have a preferred direction of propagation defined by the magnetic 
field, the effective cut-off frequency is lowered by the cosine of the inclination angle 
with respect to the vertical. This allows evanescent waves to propagate. It should be 
recalled, however, that mainly vertical propagation is observed in the photosphere 
and chromosphere in plage regions (see e.g. Krijger et al., 2001, Centeno, Collados, 
and Trujillo Bueno, 2006), thus it is unclear whether the mechanism suggested by 
De Pontieu, Erdelyi, and Stewart (2004) is at work. Alternatively, the decrease of 
the effective acoustic cut-off frequency can be produced by taking into account 
the radiative losses of acoustic oscillations (Roberts, 1983; Centeno, Collados, and 
Trujillo Bueno, 2006). If the radiative relaxation time is effectively small (as expected 
in small-scale magnetic structures, transparent to radiation), the cut-off frequency 
can be reduced and evanescent waves can propagate. 

All works cited above apply the theory of acoustic waves in idealized atmospheres 
in order to explain the observed oscillation properties in magnetic structures. In 
magnetized atmospheres, different wave types can exist and different modes can 
be observed depending on the magnetic field configuration and the height of the 
transformation layer (cs = va) relative to the height of the formation of the spectral 
lines used in observations. Bogdan et al. (2003), Rosenthal et al. (2002), Hasan and 
Ulmschneider (2004), Hasan et al. (2003, 2005), addressed this question in their 
simulations of magneto- acoustic waves in non-trivial magnetic configurations. They 
pointed out the importance of defining the mode type and the mode conversion height 
when interpreting the observations. Due to mode mixing at the transformation layer, 
no correlation may be observed between perturbations below and above the cs = va 
height. 

All these arguments point of towards the need to have realistic, self-consistent, 
and systematic models of the wave spectrum in magnetic structures. It is necessary 
to include non-trivial magnetic field configurations in order to allow the existence 
of different mode types and mode conversion. Non-linearities should be taken into 
account for waves propagating into the chromosphere. The thermal structure should 
be realistic in order to make possible the reflection of waves with frequencies below 
the cut-off. Radiative damping also plays an important role, in particular to estimate 
the amount of energy that may deposited in the higher atmosphere due to waves. In 
the present work we address several of the above questions. We perform non-linear 
numerical modeling of magneto-acoustic waves in small-scale magnetic flux tubes 
with internal structure. We consider waves with realistic periods of three and five 
minutes, which has never been done before for non-trivial magnetic configurations. 
We study the mode transformation and oscillation spectra at different heights from 
to 2000 km as a function of the type of the driver and its period. At present, radiative 
damping of oscillations is not taken into account. This question will be addressed in 
a separate paper. 



2. Numerical Model 



We solve the basic equations of the ideal MHD, written in conservative form: 

f + v-(py) = o, (1) 
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Figure 1. Second-order flux tube magnetostatic model, calculated after Pneuman, Solanki, 
and Stenflo (1986). The tube radius is equal to 100 km in the photosphere. Left panel: magnetic 
fleld lines and contours of plasma f3 (deflned as c^/v\), with values indicated on each contour. 
Right panel: distribution of the magnetic fleld strength. 



^+V-[{E + P+^)V-B{^)]=pV-g + pQ, (3) 

^ = V X (V X B) , (4) 

where I is the diagonal identity tensor and E is the total energy: 

All other symbols have their usual meaning. The numerical code solving these equa- 
tions was used previously in Khomenko and Collados (2006) and is described in 
details in that paper. For the present work, the code has been improved to solve the 
fully non-linear equations. One of the main differences between our code and any 
other standard MHD code is that we solve equations for perturbations (including all 
non-linear terms), the magnetostatic equilibrium being analytically subtracted from 
the above system of equations. 

The term describing energy losses, (Q), is set to zero in the present study. A 
Perfect Matching Layer (PML) boundary condition (Berenger, 1994) is applied to 
the sides and top boundaries of the simulation domain. It is found to perform 
rather well in our simulations. Tests have shown that wave reflection from the upper 
boundary can be reduced to a few percent even for shock waves with Mach number 
of two or three. This has allowed us to perform longer temporal series of simulations, 
not stopping them at the time when the fastest wave reaches the upper boundary. 
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Table 1. Characteristic parameters of the model flux 
tube taken at the flux tube axis. 



Parameter 


z=0 km 


z=2000 km 


Temperature (K) 


6420 


7175 


Magnetic fleld (G) 


743 


37 


Sound speed (km s"-*^) 


8 


10 


Alfven speed (km s~^) 


4 


394 




4 


7x10-4 


Pressure scale height (km) 


151 


228 


Cut-off period (s) 


230 


280 



and to analyze these simulations in a statistical way, similar to observations. The 
boundary conditions and other numerical details are described elsewhere (Khomenko 
and Collados, 2006). The simulations analyzed below are done in two dimensions. 

At the bottom (photospheric) boundary, we specify either a vertical or horizontal 
velocity as a function of time and horizontal coordinate (x): 

Vx,z{x, t) = Vo sm{ut) X exp(-(x - xof/2a^) . (6) 

The perturbation is localized in the x-axis and has a Gaussian shape in this direction. 
The horizontal size of the pulse cr = 160 km, initial amplitude at the photosphere 
Vo = 200 m s~^ and xq corresponds to the tube axis. 

3. Magnetostatic Model 

In order to construct an initial magnetostatic structure, we have followed the method 
described in Pneuman, Solanki, and Stenflo (1986). The equations for magnetostatic 
equilibrium have been written in cylindrical coordinates. All the variables were ex- 
panded in a power series in the horizontal radial coordinate r. The series have been 
truncated after the second order in r. After some manipulations with the equations, 
one single equation can be obtained for the vertical stratification of the magnetic field 
strength at the flux tube axis (for details, see Pneuman, Solanki, and Stenflo, 1986). 
This equation has been solved numerically using an iterative interpolation scheme 
from Korn and Korn (2000). Once the magnetic field at the axis is obtained, all of the 
other variables can be calculated from analytical expressions. The free parameters of 
the model are: magnetic field strength and flux tube radius at the base of photosphere 
and filling factor. The latter parameter allows us to produce a smooth merging of flux 
tubes at some chromospheric height. By varying these parameters, we have obtained 
a solution with desired properties concerning the height of the transformation layer 
cs = ^A^ inclination of field lines, and size of the tube. The temperatures inside 
and outside the flux tube were kept constant at all heights. The initial stratification 
of the thermodynamic parameters were taken from the VAL-C atmospheric model 
(Vernazza, Avrett, and Loeser, 1981). Table 1 and Figure 1 summarize the magnetic 
and thermal properties of our magnetostatic solution. 

Since our model flux tube is a second-order approximation to a thin flux tube, it 
has horizontal variations of magnetic field strength, pressure, etc. This approximation 
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also gives rise to horizontal and vertical variations of the plasma /3 (note that through 
the present paper we define f3 as c^/v\). The transition between the outside and 
inside atmospheres is rather sharp. Due to pressure balance, the gas pressure is 
smaller inside than outside the flux tube. The magnetic field lines fan out until 
a certain height, where nearby flux tubes merge together. After that height, the 
magnetic field becomes vertical and almost homogeneous in the horizontal direction. 

The magnetostatic model shown in Figure 1 was perturbed by either a vertical or 
a horizontal driver, as follows from Equation (6). The results of the simulations for 
different driving periods are described in the remainder of the paper. 

4. "Mode Conversion" Nomenclature 

There is a varying usage of the term "mode conversion" both in astrophysics and in 
plasma physics literature. 

On the one side, the term "mode conversion" is used when the mode changes from 
slow to fast or the other way around (see, e.g. Cally (2006)). The transformation of 
a slow (magnetic) mode into slow (acoustic) mode is not called "mode conversion" , 
but instead "the changing of the character of slow mode" (same for the fast mode). 

On the other side, in a later paper by Cally (2007), the terms "conversion" and 
"transmission" are used. The term "conversion" is used to describe the change from 
acoustic to magnetic nature of the wave or the other way around, i.e. slow-slow or 
fast-fast. The term "transmission" is used to describe the change from slow to fast 
or from fast to slow. This terminology used in parts of the plasma physics literature 
as well. 

Since there seems to be no standard to describe what a mode conversion is, we 
find it necessary to define the nomenclature we utilize below in the paper. We use 
the term "mode conversion" in a more general sense. When the propagation speed of 
wave changes from fast to slow, we call it "mode conversion from fast mode to slow 
mode" (the same for the opposite) . When the physics of the restoring force changes 
from magnetic to acoustic, we call it "mode conversion from slow (magnetic) mode 
to slow (acoustic) mode" (the same for the fast mode). 

5. Horizontal Driving at 50 Seconds 

Simulations of short-period {i.e. short- wavelength) waves have the advantage that 
they make possible to understand easily what types of wave modes are generated 
and what transformations take place. 

Figure 2 gives a snapshot of the transverse and longitudinal velocity components 
at an elapsed time of 340 seconds after the start of the simulations, and Figure 3 
gives the temporal evolution of this simulation during the early stages. At each point, 
the transverse component is normalized to ^/vXpo, and the longitudinal component 
is normalized to y^cspo- This normalization allows us to compare the energy fluxes 
in both components for the case of the linear, short- wavelength perturbations. 

The periodic horizontal driver located at the lower boundary in the high-/3 regime 
generates primarily a slow magnetic mode, detected in the transverse velocity per- 
turbation inside the flux tube at lower heights. Due to the Gaussian shape of the 
driver, the fast mode is also generated. It can be seen at the early snapshots in 
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Figure 2. Snapshot of the simulations excited with a horizontal driver at 50 seconds. The 
panels show the transverse (perpendicular to the magnetic field) and longitudinal (parallel 
to the magnetic field) velocity components at an elapsed time of 340 seconds after the start 
of the simulations. The transverse component is normalized to ^/vaPo and the longitudinal 
component is normalized to y^cspo in order to highlight the variations of the energy fiux as 
a function of height. The black lines represent the magnetic field lines. The two red lines are 
contours of constant c^/v\, the thick line corresponding to va = cs and the thin line to 
= 0.1. Note the change of the position of both curves, compared to Figure 1, produced 
by wave motions. 

Figure 3 (up to 240 seconds after the start of the simulation) as a weak perturbation 
in longitudinal velocity propagating rapidly up and having a longer wavelength. It 
can be seen that the fast wave is rather weak and does not play a major role at the 
stationary stage of the simulation. Its contribution to the total energy flux at the 
upper atmosphere is only minor. Note, also, that the fast wave is present outside the 
flux tube in the field-free atmosphere. This is due to the fact that the width of the 
driver (cr) is larger than the radius of the flux tube [Equation(6)]. 

The slow mode generated at the bottom propagates vertically upward chang- 
ing the curvature of the wave front due to the horizontal gradients in the Alfven 
speed inside the flux tube. Note that there is no corresponding perturbation in the 
longitudinal velocity inside the tube. This confirms that the slow mode is purely 
transverse. In addition to this slow magnetic mode, horizontal motions at the lower 
boundary also generate a surface mode. It can be seen very well in the snapshot of the 
longitudinal velocity. This mode is localized on the magnetic-non-magnetic interfaces 
and follows the shape of the tube. The perturbation associated with the surface 
mode decays exponentially away from the tube boundary inwards and outwards. 
The surface mode is longitudinal and can be distinguished in the simulations until 
the height where flux tubes merge and the magnetic field and other parameters 
become homogeneous. 

In the vicinity of the cs = va layer, both the slow magnetic mode and the surface 
mode transfer most of their energy to the slow acoustic mode. This transformation 
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Figure 3. Time series of snapshots of the transverse (right) and longitudinal (left) velocities 
in the simulation with a horizontal driver at 50 seconds, showing the early stages of the 
evolution. Numbers give the elapsed time in seconds since the start of the simulation. The 
transverse component is normalized to ^/vaPo and the longitudinal component is normalized 
to y^cspo- Each snapshot includes 900 km in horizontal and 2000 km in vertical directions. 
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Figure 4. WKB solution for the simulations with a horizontal driver at 50 seconds. The 
format of the figure is the same as Figure 2. Blue contours are the positions of constant phase 
of the WKB slow mode solution at the different equidistant time moments. 



can be appreciated from the time evolution of the simulation shown in Figure 3, where 
it is clear that significant perturbation in the longitudinal velocity in the upper layers 
inside the flux tube appears only after the slow mode reaches the transformation 
layer. Some of the energy of the slow magnetic mode also goes into the fast magnetic 
mode at the va > cs region. As shown by the sequence of snapshots, this fast 
mode is refracted back to the photosphere. However, its energy is small compared 
to that of the slow acoustic mode in the va > cs layers. The ratio of the transverse 
and longitudinal velocities at these heights is close to 0.1. This situation should 
be contrasted with the simulations of waves in sunspots (Khomenko and Collados, 
2006), where a significant part of the energy of the driver was returned back to 
the photosphere by the fast magneto-acoustic mode. The slow acoustic mode in the 
Va > Cs region is longitudinal and can be observed at heights above 1000 km in the 
longitudinal velocity snapshot. It propagates along the magnetic field lines straight 
up to the chromosphere and carries most of the energy of the driver to the upper 
heights. Motions on opposite sides of the flux tube are in anti-phase, as a natural 
consequence of the horizontal motions of the tube at the photospheric boundary. 

The flnding that the slow magnetic mode is effectively transformed into the slow 
acoustic mode is different from the conclusions of, e.g., Bogdan et al. (2003) who 
obtain from their simulation that a signiflcant part of the slow-mode energy goes 
into the fast mode, thus preserving the magnetic nature of the wave at all heights. 
In order to better understand our simulations, we present in Figure 4 the WKB 
solution for the slow mode. The WKB approach is linear and assumes that the 
wavelength of the perturbation is much smaller than the characteristic scale of the 
variations of the background atmospheric parameters, and thus can be applied in the 
case of the short-wavelength simulations with the 50-seconds driver. The equations 
describing the WKB approach in this particular case can be found in Khomenko 
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2000 



HEIGHT, KM 

Figure 5. Height dependence of phase speeds of the slow mode (blue dotted line), surface 
mode (green dotted line) and acoustic mode (red dotted line) measured from the simulations 
compared to the characteristic propagation speeds calculated in the magnetostatic model (solid 
lines). The values of the Alfven and sound speed (va and cs) are taken at the flux tube axis. 
The value of the tube speed (ct) is computed along the trajectory of the surface mode. 

and Collados (2006). The blue contours in Figure 4 describe the position of the wave 
front at different time moments separated by equal intervals. The initial wavefront 
at the base of the atmosphere is assumed planar, i.e. the horizontal wave vector kx is 
equal to zero. The WKB solution describes rather well the position of the wavefront 
at every time moment, i.e. the waves propagate at the slow-mode speed. As can be 
appreciated from the figure, as the wave propagates upwards, the horizontal gradients 
of its phase speed produce the deformation of the wave front, so that its central 
part propagates faster. It means that at greater heights the wave vector (/c) is no 
longer directed along the vertical, but forms a significant angle with respect to the 
magnetic field vector Bq. According to Cally (2006), if the angle between k and 
is close to zero, the most effective transformation is from the slow to the fast mode. 
However, as this angle increases, the slow to slow mode transformation becomes more 
important. Since in our simulations the slow mode reaches the transformation layer 
with k significantly inclined with respect to Bq, it explains the effectiveness of the 
slow to slow mode transformation. This does not contradict the fact that the wave 
energy may propagate longitudinally, which is a result of the strong anisotropy of 
the medium. 

The slow acoustic mode develops shocks with amplitudes of about 5 km s""*^ that 
propagate straight up into the chromosphere. Due to the clear non-linear behavior 
of the acoustic mode at the greater heights, secondary harmonics are generated in 
addition to the frequency of the original perturbation of the driver. However, the 
maximum power at all heights is kept at the driving frequency. 

Figure 5 gives the characteristic speeds of the waves in the simulation. Using 
the time series of the simulations, we have calculated the phase speed of the slow 
magnetic mode (at heights from to about 800 km, where this mode exists), the 
slow acoustic mode (at heights from 500 km upward) and the surface mode. In the 
latter case, the phase speed along the curved trajectory of the surface mode was 
calculated. The retrieved phase speeds are displayed in Figure 5 as a function of 
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height. The stratifications of va, cs and tube speed (c^ = c|i;^/ (c| -\-v\)) are also 
presented. 

The comparison of the dotted and solid line curves in Figure 5 leads to coherent 
results. The slow magnetic mode propagates with the Alfven speed {va)- The surface 
mode propagates with the tube phase speed (ct) up to 400 km. Despite the fact that 
the tube and the Alfven velocities are rather similar in the lower photosphere, the 
phase speed of the surface wave is distinct from that of the slow magnetic wave. The 
former propagates with its own characteristic speed, as suggested from the analytical 
theory (see, e.g. Roberts, 1981). At heights about 700 km (where va ~ c^) the speeds 
of all modes become close to one another and the energy can be easily transferred 
between the different wave types. It can be seen that the slow magnetic mode and 
the surface mode give their energy to the slow acoustic mode. As the slow acoustic 
mode propagates upwards, shock formation occurs and the propagation becomes 
supersonic. The phase speed of the shocks increases with height in agreement with 
the increase of their amplitudes. 

The following items summarize the most important features of the simulations 
with the horizontal driver at 50 seconds: 

• Slow to slow mode transformation at the va = cs height. 

• Left-right acoustic shock wave pattern in anti-phase. 

• Surface mode excitation. 

• The period of the waves is maintained at all heights. 

The simulation with this 50 seconds driver follows the spirit of the calculations by 
Rosenthal et al. (2002), Bogdan et al. (2003) and Hasan et al. (2005), where small 
scale flux tubes were perturbed by horizontal high-frequency motions at the lower 
boundary. Similarly to these works, we also find that such motions generate acous- 
tic shocks at chromospheric heights. The differences in the wave behavior between 
these works and ours are mainly due to the difference in the initial magnetostatic 
situation. Once more, this suggests that the background magnetostatic state should 
be taken into account while interpreting observations, since different wave types can 
be generated. 



6. Horizontal Driving at 180 Seconds 

Here, we show the results of the simulations using a 180-seconds driver, still above the 
cut-off frequency. For this reason, it is not surprising that we get a behavior similar to 
the previous case. The wave pattern that appears in the simulations is basically the 
same as in the case of the 50 seconds driver, but the wavelength of the perturbations 
is larger. The driver mainly excites in the lower layers a transverse slow magnetic 
mode inside the flux tube, together with a surface mode at the interface with the non- 
magnetic surroundings. As before, the surface mode is mainly longitudinal. Figure 6 
shows some snapshots of the transverse and longitudinal velocities during one period 
of the wave after the simulations reached the stationary state. In contrast to the case 
of the short- wavelength 50 seconds driver simulation, we find the normalization of 
the velocity components to -^cspo or y^vXpo to be not appropriate for the 180 
and 300 seconds driver simulations since the group velocity of waves follows a more 
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Figure 6. Time series of snapshots of the transverse (right) and longitudinal (left) velocities 
in the simulation with a horizontal driver at 180 seconds. Numbers give the elapsed time in 
seconds since the start of the simulation. After removing the average velocity increase with 
height produced by the density fall-off, all images have the same scale. Each snapshot includes 
900 km in horizontal and 2000 km in vertical direction. 



khomenko.tex; 2/02/2008; 4:30; p. 13 



E. Khomenko, M. Collados and T. Felipe 




5 10 15 20 5 10 15 20 25 30 

TIME, MIN FREQUENCY [mHz] 




TIME, MIN FREQUENCY [mHz] 

Figure 7. Examples of the temporal evolution of the longitudinal velocity (left) and their 
corresponding power spectra (right) at different heights and positions inside the flux tube for 
the simulation with horizontal driving at 180 seconds. The coordinates are indicated on the 
right panels. 



complicated behaviour than either cs or va due to the larger wavelength and non- 
linear effects. Instead we normalize both velocity components with an exponential 
function (p^^) that allows a better visualization of the results. 

Around the Cg = va layer, mode conversion takes places and most of the energy 
goes again to the slow acoustic mode. The antisymmetric pattern with respect to the 
flux tube axis is still present in longitudinal velocity above this layer, although the 
antisymmetry is not as perfect as in the 50-seconds driver case. Shocks clearly develop 
during the evolution. Very clear examples are apparent at time 379 in the left part 
of the flux tube, and at time 452 in the right part. These non-linear effects are more 
pronounced in the higher layers. The shock formation occurs above a height of about 
1000 km. The temporal evolution of the transverse velocity shows the refraction of 
the fast magnetic mode at the higher layers where va> cs- The horizontal motions 
of the flux tube field lines in the photosphere and chromosphere due to this mode 
can also be seen during the time evolution. 

In Figure 7, the temporal evolution of the longitudinal velocity at different points 
in the flux tube, corresponding to different heights and horizontal distances from 
the axis, is shown. The first two upper panels on the left correspond to two points 
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located near the boundary of the flux tube, one on the right and the other on the 
left of the axis, 200 km above the driver. The almost sinusoidal oscillation pattern 
is clear, since non-linearities have not had time to develop. When comparing the 
position in the horizontal axis of maxima and minima for the velocity time series 
from both points, it is apparent that both points oscillate in anti-phase. The two 
middle panels on the left correspond to two points located at the same distance from 
the axis than the previous ones, but at a larger height (704 km). This is the height 
where va ~ cs and the wave transformation starts. The amplitude of the oscillation 
has increased by almost a factor of four. It can be seen how the wave takes some 
time to reach this layer. The perturbation reaches this layer about four minutes after 
the driver starts. Up to there, non-linearities are still negligible, and the sinusoidal 
aspect of the wave is still kept, as well as this antisymmetry of the velocity pattern. 
At higher layers (1664 km in the plots) velocity discontinuities appear, indicating 
the development of shocks, with amplitudes of several km s""*^ (two bottom panels 
of Figure 7). The corresponding power spectra of this velocity temporal series are 
plotted in the right-hand part, and they all demonstrate that the frequency of the 
driver (5.56 Hz) remains unaltered at all heights. The presence of several harmonics 
in higher layers is indicative of the non-linearities giving rise to the saw-tooth velocity 
profiles, rather than sinusoidal. Note that the shocks are in antiphase at the both 
sides of the tube. 

Probably, the most important difference with respect to the previous case with 
a driver of 50 seconds, lies in the fact that the slow-to-slow mode conversion is 
not as efficient now (see Cally, 2006), and a larger part of the energy of the slow 
mode excited by the driver is transferred to the fast (transverse, and hence, mainly 
magnetic) mode. The atmosphere above the va = cs oscillates uniformly to left and 
right with small phase differences. The horizontal motions are in phase in the whole 
atmosphere due to the large oscillation wavelength. 

Summarizing, the most important features of the simulations with the horizontal 
driver at 180 seconds are: 

• Partial slow to slow (acoustic) transformation at va = cs height. 

• Horizontal flux tube motion at 2000 km due to the fast (magnetic) wave. 

• Left-right shock wave pattern in anti-phase. 

• Surface mode excitation. 

• Period of the wave maintained with height. 

7. Horizontal Driving at 300 Seconds 

Figure 8 shows the temporal evolution of the transverse and longitudinal velocities 
over one wave period, after the simulations reach the stationary state, and Figure 9 
gives the temporal variation of the longitudinal velocity at selected points, together 
with their corresponding power spectra, for the simulations with a horizontal driver 
with a period of 300 seconds. The frequency of the driver is below the cut-off fre- 
quency. As in two previous cases, the driver excites the slow and surface modes. 
Note that, unlike the shorter-period simulations, there are significant motions ex- 
cited outside the flux tube in the non-magnetic atmosphere as well. The nature of 
these motions is partly acoustic and partly convective since the lower layers of the 
photosphere are still unstable to convection. 



khomenko.tex; 2/02/2008; 4:30; p. 15 



E. Khomenko, M. Collados and T. Felipe 




Figure 8. Time series of snapshots of the transverse (right) and longitudinal (left) velocities 
in the simulation with a horizontal driver at 300 seconds. Numbers give the time elapsed in 
seconds since the start of the simulation. After removing the average velocity increase with 
height produced by the density fall-off, all images have the same scale. Each snapshot includes 
900 km in horizontal and 2000 km in vertical direction. 
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Figure 9. Examples of the temporal evolution of the longitudinal velocity (left) and their 
corresponding power spectra (right) at different heights and positions inside the flux tube for 
the simulations with horizontal driving at 300 seconds. The coordinates are indicated on the 
right panels. 



The slow-to-slow mode transformation is only partial, and part of the energy of 
the slow mode is transferred to the fast mode above the va = cs layer, producing a 
horizontal shaking of the tube that propagates upwards with the Alfven speed. The 
behavior of the fast mode seen in the transverse velocity above va = cs height is 
similar to the three-minute period simulations, but the left-right symmetry pattern 
is lost. 

Most of the energy in high layers is nonetheless still carried by the slow (acoustic) 
mode, giving rise to shock waves with amplitudes of 10 — 15 km s~^ above 1000 km 
height. The antisymmetric nature of these shocks has disappeared. Left and right 
parts of the tube are almost in phase now. However, the most conspicuous difference 
lies in the period of the wave. As it corresponds to an evanescent wave, the five-minute 
oscillation is rapidly damped, and the residuals coming from non-linear effects at 6 
Hz, twice the driver frequency, already dominate at a height of 700 km. This double 
frequency is clear in the temporal series of the velocity, and in the corresponding 
power spectra, at this height (two middle panels in Figure 9). From 700 km upward, 
the velocity increases again, due to the density decrease, giving rise to shocks in the 
upper part of the simulation domain. The main frequency is maintained around 6 
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mHz, once the five-minute oscillation has been damped. Kalkofen et ai (1994) and 
Fleck and Schmitz (1991) showed that the change of the dominant period of acoustic 
oscillations with height, from five to three minutes, is due to the resonant excitation of 
waves at the atmospheric cut-off frequency {i.e. the cut-off frequency corresponding 
to the temperature minimum). However, this effect may not be dominant in our 
simulations since the acoustic mode only appears above 1000 km, i.e. well above 
the temperature minimum. At heights above 1000 km, the cut-off frequency is again 
around three mHz and the resonant excitation cannot produce oscillations at the 
higher frequency. We conclude that the generation of oscillations at the first harmonic 
of the driver in our simulations is mostly a non-linear phenomenon. 

As in the previous sections, the most important properties of the wave propagation 
from the simulations with the horizontal driver at 300 seconds are: 

• Partial slow to slow (acoustic) transformation at va = cs height. 

• Horizontal flux tube motion at 2000 km due to the fast (magnetic) wave. 

• Left-right shock wave pattern in phase. 

• Surface mode excitation. 

• Period of the wave not maintained with height. 



8. Vertical Driving at 300 Seconds 

Finally, the results after exciting the flux tube with a 300-seconds vertical driver are 
presented. Figure 10 gives some snapshots of the velocity evolution during one period 
with the same format as in the previous figures. The time series shown in the figure 
is taken after the simulations reach the stationary state. Now, the fast (acoustic) 
wave is excited directly by the driver inside the flux tube in the deeper layers (see 
lower panels of Figure 10). The surface wave exists as well at the magnetic/field-free 
interface. There are some acoustic disturbances that can be observed in the field-free 
atmosphere surrounding the tube at the lower layers. After some time has elapsed 
since the start of the simulation, convective vortices appear in the non-magnetic 
atmosphere. 

At the Va = cs layer, the fast to slow transformation occurs, keeping the acoustic 
nature of the wave. Almost no energy is transferred to the fast magnetic mode in high 
layers. It can be seen in the upper panels of Figure 10 that the transverse velocity 
inside the flux tube remains close to zero all of the time. The magnetic field lines 
are compressed and expanded as the wave passes through the atmosphere, but no 
horizontal motions are observed. 

The perturbation is almost perfectly symmetric at all heights with respect to 
the tube axis, as demonstrated in the four upper panels of Figure 11, where the 
time series in symmetric positions inside the flux tubes are given. The most striking 
fact comes from the different behavior of the velocity at points near the tube axis 
or at the tube boundary in the upper heights. At deep layers, the propagation is 
mainly linear showing the five- minute sinusoidal pattern of the driver (see the 
upper panels in Figure 11). From there upwards, two shock waves develop at the 
beginning of the series, which are followed by a highly damped oscillation. The 
period of these shocks is three minutes. After the simulation reaches the stationary 
state, the dominant periods are different at the tube axis and off axis. For points near 
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Figure 10. Time series of snapshots of the transverse (right) and longitudinal (left) velocities 
in the simulation with a vertical driver at 300 seconds. Numbers give time in seconds since the 
start of the simulation. After removing the average velocity increase with height produced by 
the density fall-off, all images have the same scale. Each snapshot includes 900 km in horizontal 
and 2000 km in vertical direction. 
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Figure 11. Examples of the time evolution of the longitudinal velocity (left) and their cor- 
responding power spectra (right) at different heights and positions inside the flux tube in 
simulations with vertical driving at 300 seconds. The coordinates are indicated on the right 
panels. The blue dotted line at the right panels gives the power spectra only for the last half 
of the simulations after reaching the stationary state. 



the axis, the five-minute oscillation does not suffer such a large amplitude damping 
and is the dominate period once the oscillation reaches a stationary state (see the 
power spectra corresponding to the second half of the simulation as a blue line in 
Figure 11). For off-axis points, there is a clear change with height of the dominating 
frequency. The contribution of the five-minute fundamental frequency is reduced as 
the wave propagates, increasing the relative importance of the first harmonic at 6 
mHz (lower panel of Figure 11). Thus, in this simulation, 3 and 6 mHz oscillations 
coexist in different parts of the tube. In any case, the amplitudes are rather small 
when compared to those derived from horizontal drivers and only weak eventual 
shocks are observed with amplitudes of 2 — 3 km s~^. 

To summarize, the most important features of the simulations with a horizontal 
driver of 300 seconds are: 

• Fast to slow transformation at va = cs height. 

• Weak eventual shocks 

• Period on- axis: five minutes 
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• Period off- axis: three minutes 

If the simulations with horizontal drivers can be understood in terms of mode 
transformation, and wave amplification/damping, this is not the case for this sim- 
ulation with a 300-seconds vertical driver. Why does the main period of the wave 
above the va = cs layer depend on the distance to the axis? Why do shocks develop 
at the beginning of the series and not later? More simulations are needed in different 
flux tubes, with varying initial conditions, to understand this complex behavior. 



9. Discussion and Conclusions 

We have performed simulations of waves in flux tubes with non-trivial magnetic 
field configuration. The model flux tube has horizontal and vertical variations of 
magnetic field strength and gas pressure and of the ratio between the characteristic 
wave speeds va and cs- These properties have allowed us to study wave propagation 
and mode transformation inside the flux tube. Waves are excited by a photospheric 
driver with a period of three and five minutes, for the first time with a magnetic 
field configuration of this kind. Different wave patterns are observed depending on 
the period and on the type of the driver. The following items summarize our findings: 

• Horizontal motions of the flux tube at the bottom photospheric boundary 
generate a slow magneto-acoustic mode inside the tube and a surface mode 
at the magnetic/field free interface. 

• After the slow magneto-acoustic mode and the surface mode reach the height 
where va ~ c^, their energy is effectively transformed into a slow acoustic 
mode in the high atmosphere where va > cs- Only a small part of the driver 
energy is returned back to the photosphere by the fast magneto-acoustic mode. 

• The slow acoustic mode propagates vertically along the magnetic field lines 
in the atmosphere where va > cs^ forms shock waves above 1000 km with 
amplitudes of 5 — 15 km s""*^ and remains always within the same flux tube. 
Thus, it can deposit effectively most part of the energy of the driver into the 
chromosphere. 

• If the frequency of the horizontal driver is above the acoustic cut-off, non-linear 
wave propagation in the tube occurs with the dominant period of the driver at 
all heights. Non-linear effects produce higher harmonics, but their amplitude is 
small. 

• If the frequency of the horizontal driver is below the acoustic cut-off, the dom- 
inant period of the longitudinal velocity changes with height from five to three 
minutes due to non-linear generation of the higher harmonics. 

• The five-minute vertical perturbation of the flux tube at the bottom pho- 
tospheric boundary generates a fast acoustic mode that propagates upward 
through the transformation va = cs layer, without changing its nature and 
forming only eventual shocks in the chromosphere. 

• After reaching the stationary state in the simulations with a vertical driver at 
300 seconds, both three and five minute oscillations co-exist inside the flux tube 
at chromospheric heights. The dominant period at the axis is five minutes and 
the off- axis dominant period is three minutes. 
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Our simulations suggest once again that the properties of waves observed in 
magnetic structures are the direct consequence of their magnetic configuration, tem- 
perature, magnetic field strength, and the height where the mode transformation 
occurs. When analyzing observations, it is important to know whether they corre- 
spond to the level below or above the layer va = cs- If the motions that excite 
oscillations inside flux tubes are purely horizontal, no correlation may be observed 
between velocity variations measured in the photosphere and the chromosphere. This 
is because only an acoustic mode can be detected in observations (at least at disc 
center) since it produces significant vertical velocities, unlike the slow magnetic mode 
existing in the photosphere. The acoustic mode is only generated above 700 — 1000 
km in our simulations with a horizontal driver. In contrast, if the driver that excites 
oscillations has a vertical component, the acoustic mode is generated already in 
the photosphere. In this case, the oscillations in vertical velocity are coherent at 
photospheric and chromospheric heights. 

In general, both simulations with horizontal and vertical drivers with 300-second 
periodicity show the change of the dominant wave period with height. This property 
is different from observations of plage and network regions where the five minute 
periodicity is preserved also in the chromosphere (Lites, Rutten, and Kalkofen, 1993, 
Krijger et a/., 2001, De Pontieu, Erdelyi, and de Wijn, 2003, Centeno, Collados, and 
Trujillo Bueno, 2006, Bloomfield et a/., 2006, Vecchio et a/., 2007). While our study 
may help to identify the wave modes observed in small-scale magnetic structures, 
further work is needed in order to explain the behavior of the oscillatory spectra 
with height. 
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